    #https://www.cnblogs.com/lc1217/p/6514734.html
##最小二乘法
import numpy as np   ##科学计算库 
import scipy as sp   ##在numpy基础上实现的部分算法库
import matplotlib.pyplot as plt  ##绘图库
from scipy.optimize import leastsq  ##引入最小二乘法算法
from chinese_calendar import is_holiday
# import matplotlib 
import datetime

today = datetime.date.today()


def n_days_before(st_day,n):
    """
    st_day:开始日期
    :type st_day: datetime.date | datetime.datetime
    """
    while(n>0):
        if not is_holiday(st_day):#if workday ,mean if open stock
            n=n-1
        st_day=st_day-datetime.timedelta(1)
    return st_day

days=7
end_datetime=datetime.date.today()
start_datetime=n_days_before(end_datetime,days) 
end_date_str=end_datetime.strftime("%Y%m%d") 
start_date_str=start_datetime.strftime("%Y%m%d") 

# start_date_str=(today-datetime.timedelta(8)).strftime('%Y%m%d') #9日交易日(含休市)数据拟合


from matplotlib.pylab import mpl  
mpl.rcParams['font.sans-serif'] = ['FangSong'] # 指定默认字体  
mpl.rcParams['axes.unicode_minus'] = False #  解决保存图像是负号'-'显示为方块的问题 

''' 
设置样本数据,真实数据需要在这里处理
'''
##样本数据(Xo,Yh),需要转换成数组(列表)形式
# Xo=np.array([6.19,2.51,7.29,7.01,5.7,2.66,3.98,2.5,9.1,4.2])
# Yh=np.array([5.25,2.83,6.41,6.71,5.1,4.23,5.05,1.98,10.5,6.3])
import tushare as ts
ts.set_token("6667cd4a2326f2f937062a0f4fb59aea5c56d13b1f6f26225f115fe9")
pro = ts.pro_api()

symbols=['002108.SZ','002600.SZ','002146.SZ','002938.SZ','002269.SZ'] #[沧州明珠,领益制造,荣盛发展,鹏鼎控股,美邦服饰]
df7 = pro.daily(ts_code=symbols[0], start_date=start_date_str)
df7["wave"]=df7["high"].values-df7["low"].values
# plt.plot(df7["wave"])
# plt.show() 
Yh=df7["high"].values
Yl=df7["low"].values
Xo=df7["open"].values
Xc=df7["pre_close"].values
# print("--------------------------")
# print(Xo)
# print(Xc)
# print("--------------------------")
date= df7["trade_date"]

'''
设定拟合函数和偏差函数
函数的形状确定过程:
1.先画样本图像
2.根据样本图像大致形状确定函数形式(直线、抛物线、正弦余弦等)
'''

##需要拟合的函数func :指定函数的形状y=kx+b
def func(p,x):
    m,k,b=p
    m=0 #此时为1次式
    return k*x+b+m*x*x

##偏差函数:x,y都是列表:这里的x,y更上面的Xo,Yh中是一一对应的
def error(p,x,y):
    return func(p,x)-y



'''
主要部分:附带部分说明
1.leastsq函数的返回值tuple,第一个元素是求解结果,第二个是求解的代价值(个人理解)
2.官网的原话（第二个值）:Value of the cost function at the solution
3.实例:Para_high=>(array([ 0.61349535,  1.79409255]), 3)
4.返回值元组中第一个值的数量跟需要求解的参数的数量一致
'''

#k,b的初始值,可以任意设定,经过几次试验,发现p0的值会影响cost的值:Para_high[1]
p0=[0,1,2]

#把error函数中除了p0以外的参数打包到args中(使用要求)
Para_high=leastsq(error,p0,args=(Xo,Yh))
Para_low=leastsq(error,p0,args=(Xo,Yl))
#https://baike.baidu.com/item/黑塞矩阵/2248782?fr=aladdin
#读取结果
m,k,b=Para_high[0]
ml,kl,bl=Para_low[0]
err=np.sum(error(Para_high[0],Xo,Yh)**2)/len(Xo)
print("平均误差",err)
print("预测目标数据:",func(Para_high[0],Xo))
print("真实目标数据:",Yh)
print("误差:",error(Para_high[0],Xo,Yh))
print("m=",m,"k=",k,"b=",b)
print("cost:"+str(Para_high[1]))
print("求解的拟合直线为:")
print("y_high="+str(round(m,2))+"x*x+"+str(round(k,2))+"x+"+str(round(b,2)))
print("y_low="+str(round(ml,2))+"x*x+"+str(round(kl,2))+"x+"+str(round(bl,2)))
xo1=float(input("请输入当日开盘价")) 
print("预估最高价:",func(Para_high[0],xo1))
print("预估最低价:",func(Para_low[0],xo1))
'''
   绘图,看拟合效果.
   matplotlib默认不支持中文,label设置中文的话需要另行设置
   如果报错,改成英文就可以
'''

#画样本点
plt.figure(figsize=(8,6)) ##指定图像比例: 8:6
plt.scatter(Xo,Yh,color="green",label="样本数据",linewidth=1) 



#画拟合直线
x=np.linspace(0,int(np.max(Xo))+1,100) ##在0-15直接画100个连续点
y=k*x+b ##函数式
plt.plot(x,y,color="red",label="拟合直线",linewidth=1) 
plt.legend(loc='lower right') #绘制图例
plt.xticks=date
# plt.show()


import statsmodels.api as sm

def 一元线性回归():
    x1 = Xo
    y = Yh

    x_=sm.add_constant(x1) #列向量前增加全是1的向量,变成n*2矩阵
    model=sm.OLS(y,x_) #Ordinary Least Squares 常用最小二乘法
    results=model.fit()

    print(results.params)

# 一元线性回归()


def 多元线性回归():
    x1 = Xo
    x2 = Xc
    y = Yh
    x_=np.column_stack((x1,x2))#注意加()参数只有1个

    x_=sm.add_constant(x_) #列向量前增加全是1的向量,变成n*2矩阵
    model=sm.OLS(y,x_) #Ordinary Least Squares 常用最小二乘法
    results=model.fit()
    d,i,j=results.params
    y_=d+i*x1+j*x2
    error=y_-y
    print("真实值:",y,"\n估计值:",y_,"\n误差:",error,"\n误差比例:",error/y)
    print(results.params)

多元线性回归()